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DcteniiininR  a  target  aircraft's  position,  velocity,  and  ac¬ 
celeration  from  measurement  of  its  position  is  a  major  problem  in  ^ 
modem  antiaircraft  gun,  laser,  and  missile  fire  control*  ffteferen CITS" 

1  and  2)V  Many  versions  of  the  optimal  or  Kalman  state  estimator  have 
been  used  with  varyinj;  degrees  of  succcs^{Refcrcnccs  3,  4,  and  5). 

'^^This  paper  offers  a  new  technique  which  allows  detailed  aircraft  flight 
d)'naraics  and  constraints  to  be  incorporated  into  the  general  Kalman 
structure. 


Tlib'^ilnian  estimator  is  viewed  as  a  set  of  algorithms  that 
uses  sensor '(d.g. ,  radar)  data,  the  statis.ical  properties  of  the 
sensor  errors,  the  equations  of  the  target  motion,  and  the  statistical 
properties  of  the  errors  of  these  equations  to  produce  estimates  of 
target  position,  velocity,  and  acceleration. 

Let  X  denote  the  target's  state  9  vector.  The  discretized 
dynamical  equation  is: 


X(n*l)  ■  fCn)  X(n)  ♦  N(n), 


(1.0) 


where  n  denote:;  the  time  step.  The  state  transition  or  plant  9X9 
matrix,  <^(n),  eKil)odles  target  d)'namics.  Ihe  plant  noise  random  9 
vector,  N(n),  accounts  for  liuth  the  aircraft  dynamics  neglected  by 
the  pKant  and  the  randomness  of  pilot  com.inds. 


/f 


jOj  if/' 


SCMCDtR 


Let  li  denote  expected  value  and  superscript  T  denote  trans- 
t  ll  ^ 

pose,  llie  n  tarp.et  state  estimate,  X(n),  and  the  corresponding 

A  A 

State  error  covariance  matrix,  P(n)  =  L(X(n)  -  X(n))(X(n)  -  X(n))  , 
are  extrapolated  to  time  step  n+I  by  means  of  equation  1.0.  The 
extrapolations  are  corrected  by  the  sensor  measurement  3  vector 
Z(n+1)  received  at  time  step  n+1. 

Hie  extrapolation  stage  consists  of  equations  1.1  and  1.2 

below, 

X®(n)  =  4  00  XOO  ♦  N(n)  (1.1) 

where  superscript  e  denotes  extrapolated  value  and  (  )  denotes  most 
probable  value.  Tlie  extrapolated  state  error  covariance  P^^Cn),  is 
defined  as, 

E(X(n+l)  -  X®(n))(X(n*l)  -  X®(n))^. 
lliis  definition,  equation  1.1,  and  the  assumption 

E(X(n)  -  X(n))  (N(n)  -  N(n)^  =  0  yields: 

P®(n)  =  <-(n)  P(n)  ♦  Q(n).  (1.2) 

where  the  plant  noise  covariance  9X9  matrix  Q(n),  is  defined  as: 

Q(n)  «=  E(Nfn)  -  N(n))(N(n)  -  N(n))’’’. 

At  time  step  n+1,  a  new  sensor  measurement  is  received, 

Z(n+1)  =  ll(n+l)  X(n+1)  +  V(n+1). 

The  state  to  measurement  transformation  3X9  matrix,  H(n+1),  relates 
the  state  vector  to  the  measurement  vector  (measurement-state  combina¬ 
tions  arc  discussed  in  Reference  6).  Tlic  zero  mean  random  sensor 
error  is  V(n*l)  (e.g.,  thermal  noise). 

Tlie  estimator  weighs  the  confidence  in  the  measurement 

Z(n+1)  against  the  confidence  in  the  extrapolated  state  X  (n)  and  dis¬ 
tributes  the  3  vector  difference,  Z(n+1)  -  H(n+1)  X(n),  to  correct 

appropriately  the  extrapolated  state  9  vector,  X  (n).  Consider, 


X 
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J  =  (X(n*l)  -  X®(n))  (X(n+1)  -  X®(n))'^ 

+  (Z(n+1)  -  ll(n+l)  X(n+1))  R"^(n*l)  (Z(n+1) 

-  Il(n+1)  X(n+1))^ 

where  R(n+1),  the  measurement  noise  covariance  matrix,  equals 
E  V(n+1)  V(n+1)^. 

A 

Choosing  X(n+1)  to  minimize  J  has  the  desirable  effect  of  forcing 

X(n+1)  close  to  X®(n)  when  P^(n)  is  small  (i.e.,  high  confidence  in 

X*^(n))  and  forcing  ll(n+l)  X(n+1)  close  to  Z  when  R(n+1)  is  small 
(i.e.,  high  confidence  in  Z).  Setting  the  derivative  of  J  to  zero 
yields: 

X(n+1)  =  X°(n)  +  K(n)  (Z(n+1)  -  H(n+1)  X®(n))  (1.3) 

where  the  Kalman  gain  9X3  matrix  K(n)  is  given  by: 

K(n)  =  P®(n)  M’^(n+1)  (H(n+1)  P®(n)  H’’’(n+1) 

+  R(n+1))”^  (1.4) 

The  state  error  covariance  matrix  at  the  n+1  time  step, 

P(n+1)  =  E(X(n-H)  -  X(n+1))  (X(n+1)  -  X(n+1))^, 

e  T 

computed  with  the  assumption  E(X(n+l)  -  X  (n))V  (n)  =  0  yields: 

P(n+1)  =  P®(n)  -  K(n)  H(n+1)  P®(n).  (1.5) 

Equations  1.1  through  1.5,  recognized  as  the  Kalman  algorithms 
(Reference  7),  may  be  used  iteratively  to  compute  X  for  all  n. 

2.  DEVELOPMENT  OF  $,  Q  AND  N 

Once  the  algorithms  have  been  initialized,  4',  R,  Q,  M,  and 
N  can  be  updated  to  incorporate  information  generated  from  past  state 
estimates.  Nonlinear  processes  always  require  this  sort  of  self- 
adapting  mechanism.  Adaptive  bandwidth,  hypothesis  testing,  and  resid¬ 
ual  testing  filters  arc  further  examples  that  use  past  state  estimates 
(Reference  8). 

"Dvo  self-adapting  methods  which  incorporate  detailed  target 
dynamics  and  constraints  into  4’,  Q,  and  N  are  developed  below. 
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2.]  The  Plant  Noise  Moan  Method. 


Express  the  tarpet  state  in  earth  coordinates  and  choose  a 
plant  matrix,  4*,  that  represents  the  dynamics  of  an  aircraft  whose 
Cartesian  earth  frame  acceleration  is  constant  over  a  time  step. 


^3x3 

dt  I_  , 
3x3 

I 

2  3x3 

♦  = 

^3x3 

^3x3 

dt  I_  . 

3x3 

®3x3 

Sx3 

^3x3 

where  dt  is  the  time  step  length,  and  I  is  the  identity  3x3  matrix. 
The  plant  noise  vector  then  takes  the  form, 


N 


“''a 

dtA, 


,  N  = 


dt^^ 
“6  ^ 

2_  3 
dtA_ 


Q  = 


dt" 

36 

c 

dt' 

12 

6 


dt' 

12 

4 

1 

2 


6 

dt' 

2 

dt' 


Q’ 


3x3 


where  A  is  the  target's  earth  frame  acceleration  rate  3  vector  and  Q' 
is  the  target's  earth  frame  acceleration  rate  covariance  3x3  matrix. 
The  above  formulation  is  essentially  a  Taylor  expansion  of  the  target' 
position  vector.  The  plant,  <5>,  accounts  for  the  zero,  first,  and  sec¬ 
ond  derivatives  while  the  plant  noise,  N,  is  the  remainder  term.  N  is 
thus  identified  with  the  target's  acceleration  rate. 


Tlic  acceleration  rate  statistics  are  controlled  both  by  the 
desires  of  the  pilot  and  the  maneuver  capabilities  of  the  aircraft. 
For  example,  the  greater  the  target's  capability  to  change  the  number 
of  g's  that  it's  pulling,  the  poorer  the  constant  acceleration 
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assumption  becomes  (i.c.,  the  poorer  the  plant  t’j  and  the  larger  the 
acceleration  rate  variance.  The  following  analysis  translates  physical 
projicrties  in  aircraft  body  coordinates  to  mathematical  pro|)crtics  in 
the  earth  reference  frame.  For  clarity,  a  zero  gravitational  field  is 
considered  first  and  the  results  then  extended  to  the  more  realistic 
Ig  field. 


At  any  particular  instant,  a  circle  can  be  fitted  to  the 
target's  flight  path  so  tliat  the  target  velocity  vector  is  tangent  to 
its  circumference  and  that  the  component  of  target  acceleration  per¬ 
pendicular  to  the  velocity  vector  points  towards  its  center  (Reference 
9).  Tlie  target  acceleration  is  then  expressed  as, 

A  =  Vly  *  kV^l„  (2.1) 

where  V  equals  the  magnitude  of  the  target  velocity,  1^  is  the  unit 
vector  in  the  direction  of  target  velocity,  and  1^^  is  the  principal 

normal.  1„,  which  is  perpendicular  to  the  target  velocity  vector  and 

^  -I 

points  toward  the  center  of  curvature,  equals  (kV)  The  curvature, 

k,  is  given  by  |Vly  cross  a|V'^  (1/k  =  radius  curvature). 

In  this  coordinate  frame,  a  simple  relationship  exists  be¬ 
tween  target  acceleration  and  physical  aircraft  quantities.  Tlie  force 
on  the  target  can  be  expressed  as, 

F  =  Tly  +  APlj^  (2.2) 

where:  F  =  force  per  unit  mass  vector, 

T  =  net  thrust  per  unit  mass,  and 
AP  =  air  pressure  force  per  unit  mass. 

Comparing  equations  2.1  and  2.2, 

T  =  V  „  ’  (2.3) 

AP  =  kV"^ 

A  simple  relationship  also  exists  between  acceleration  rate  and  the 
pilot  controlled  physical  quantities  of  thrust  rate,  roll  rate,  and  g 
rate.  Differentiating  equation  2.1, 

A  =  fly  .  Ti^,  *  APl^  *  APijj  (2.4) 

and  noting  the  Frenet  equations  (Reference  9) , 


sciu;ui:r 


«=  -kVl  ♦  tVI 
N  V  B 


(2.5) 


whore  x  is  the  torsion,  iV  =  K  =  target  roll  rate,  and  cross 

Ijg,  and  resolving  A  along  the  F'rcnet  axis,  yield: 


•  *23 

A  dot  ly  =  T  -  k  V 

A  dot  L,  =  kVV  ♦  a‘p 
A  dot  1„  =  kV^W 

D 


(2.6) 


t,  Al’,  and  W  arc  pilot  controlled  and  are  modeled  as  zero  mean  random 
variables.  Thcrcfoi'e,  in  the  reference  earth  coordinate  frame. 


A  = 


(2.7a) 


T  O'  T 

FR  ^  F  FR 

/V:  0.0 

0.0 

(  0.0  n  aV^ 

0.0 

\0.0  0.0 

kV^EK 

(2.7b) 


where  A  is  recognized  as  the  most  probable  value  of  the  acceleration 
rate  in  reference  coordinates,  and  Q'  as  the  target's  earth  frame 
acceleration  rate  covariance.  Tpj^  is  the  Frenet  frame  to  reference 

frame  coordiantc  transformation  and  Q'  the  acceleration  rate  covari- 

r 

ance  in  Frenet  coordinates.  The  random  variables  T,  AP,  and  W  are 
assumed  independent. 

N,  N,  and  Q  arc  computed  in  terms  of  V,  k,  Tj^p,  £(7*^),  E(AP^), 
and  E  .  V,  k,  and  arc  computed  each  time  step  in  terms  of  the 

Kr 

most  recent  state  estimate  (T^p  =  Tpj^) , 
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-  2  *  2  “  2  ^ 

V  =  (X/  4 

V  «  (X^  X^  ♦  Xj  Xg  4  x^  Xg)/V 

k  =  ((Xg  Xg  -  Xg)2  4  (X^  Xg  -  x^ 
4  (x^  Xg  -  Xg 


Trf“  =  V'^ 

T  =  X  /V 

IrP  Ag/V 

T  =  X  /V 

^RF  ^ 


T„/'  =  (X,  -  V/“n 

Trf''  =  (*8  -  *5  \’/‘=N 

V’  =  t*9  -  *6  N’/'n 


31 

.P  12 

.P  23 

13 

^  22 

“  ^RF 

T 

RF 

.  T 

‘rf 

^RF 

32 

a.  21 

.P  13 

a.  11 

.P  23 

“  ^RF 

T 

‘rf 

■  ^RF 

^RF 

33 

11 

a-  22 

a.  21 

n.  12 

“  ^RF 

^RF 

■  ^RF 

^RF 

^RF  "  %  ^RF  ■  'rF  ^RF  '•  ■' 

AAAAAA^  2 

where  A.-  =  (X.  X_  4  X-  X„  4  X,  Xo)/V  and  G.,  =  kV^.  Typical  maximum 
1  4/boDy  N 

aircraft  rates  are  given  in  Table  2-1.  Values  used  during  validation 

of  the  method  are  also  given. 

TABLE  2.1  AIRCRAFT  ACCELERATION  RATE 


Maximum  Rates 


T  (M/sec^)  10.0 
AP  (M/sec^)  20.0 
W  (rad/sec)  3.2 


Mechanized  Values 
(One  Standard  Deviation) 


Extending  this  idea  to  a  gravitational  field  requires  distri¬ 
bution  of  gl^  to  the  Frenet  axes. 
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where: 


T  =  V  + 

AP  =  ^  .  C^-9) 

T  =  target  thrust  minus  drag,  and 
AP  =  aerodynamic  force  perpendicular  to  the  target's 
wings,  and  -gl^  is  the  grativational  acceleration 
vector. 


Equation  2.9  implies; 


where : 


V  =  T  -  gkVl^dotlj^ 

(kV^)  =  CjAP  +  C2W  +  gkVl^dotly 
Cj  =  AP/fgl^dotlj^  +  kV^) 

C2  =  -gkV^l^dotlg/Cgl^dotljj  +  kV^). 


Equation  2.6  becomes: 

Adotly  =  t  -gkVl^dotlj^  -  kV 
Adotl„  =  kVV  +  gkVl_dotl,,  +  c.AP  +  c.W 

N  L  V  1  ^ 

Adotl„  =  kV^W 

D 


Hence,  analogus  to  equation  2.7: 


/  -gkV  l^dotl^  -  k^V^  \ 
A  =  \  kVV  +  gkV  l^dotl^  ) 


0.0 


ET^  0.0 

0.0 

2*2  22 

,  .,2 

0.0  Cj  EAP  +  C2  EW 

C2l^V^EW 

,  ..2  2 

.24  2 

0.0  C2kV  EW 

k  V  EW 

(2.10) 
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2.2  llic  Dynamical  Pl;mt  Method. 

The  acceleration  rate  mean  can  readily  be  incorporated  into 
the  plant  matrix  '}'.  Let  Aj.  equal  the  acceleration  vector  in  the 
Frenet  frame,  then, 


\  “  ^FR  ^F  *  *FR  ^F 


which  can  be  written; 


where : 


\  ^FR  ^^FR  ^*FR^  ^FR  \  ^  ^FR  ^F 


T‘ 

T  T 
‘fR  FR 


Therefore, 


\  “  ^FR 


FR 


+  T 


FR 


AP 

-kV^W' 


(2.11) 


is  the  desired  zero  mean  plant  noise  formulation.  For  example,  con¬ 
sider  an  aircraft  flying  in  a  deterministic  circle,  then,  .its  X  and 
Y  components  of  acceleration  satisfy:, 

=  G  cos  CO  t 

Ay  =  G  sin  (0  t 

Aj^  =  -0)  G  sin  (0  t 

Ay  =  (0  G  cos  (0  t  (2.12) 

where  R  =  lA  is  the  radius  of  the  flight  path  circle,  lo  =  V/R,  and 
G  =  |a|.  Equation  2.12  can  be  written  as  a  special  case  of  equation 
2.11. 
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The  example  also  illustrates  the  essential  equivalence  of  the 

two  methods.  Since  T  =  aV  =  W  -  V  =  0  and  =  Gto,  equations  2.7a 

and  2.8  imply, 


=  T 


TR 


Therefore, 


A  = 


which  is  consistent  with  equation  2.12. 


Both  methods  essentially  color  the  target's  acceleration  by 
fitting  the  instantaneous  Frenet  circle  to  the  latest  estimated  target 
state.  With  proper  account  given  to  the  gravitational  field,  this  ex¬ 
presses  the  plant  noise  in  the  target  related  quantities  t,  a’p,  and  W. 


3.  CONCLUSIONS 


Through  the  use  of  the  Frenet  equations,  target  roll  rate, 
thrust  rate,  and  g  rate  arc  embedded  into  the  estimator  formulation. 
This  enables  the  estimator  to  provide  more  accurate  state  error  co¬ 
variances  resulting  in  improved  state  estimates.  Similarly,  incor¬ 
porating  the  deterministic  portion  of  the  target's  acceleration  rate, 

N,  provides  a  further  improvement,  particularly  helpful  along  coordin¬ 
ated  target  turns. 

Computer  simulations  were  used  to  compare  the  plant  noise 
mean  estimator  to  a  close  copy  of  the  Gun  Low  Altitude  Air  Defense 
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System  (CII.MUS)  estimator  (Reference  4).  The  GLAADS  estimator  was 
chosen  for  reference  because  it  performed  very  successfully  during  the 
recent  f'eld  tests.  It  takes  excellent  advantage  of  modem  analytical 
techniques  and  most  imijortantly  it  is  formulated  similar  to  the  plant 
noise  mean  estimator.  'I'he  principal  difference  is  that  the  GLAADS 
estimator  lacks  an  N  term  and  uses  a  simpler  plant  error  covariance. 

The  differences  in  pcrfonn.-ince  arc  quantified  in  Figure  3.1  which  com¬ 
pares  the  magnitude  of  the  acceleration  errors  of  the  two  estimators 
for  an  aircraft  in  a  4g  coordinated  turn.  After  a  2-3  second  settling 
time,  the  advantage  of  including  detailed  target  dynamics  is  apparent. 

Additional  benefits  arise  naturally  from  the  physical  approach 
taken  by  the  plant  noise  mean  and  dynamical  plant  methods.  For  example. 
Table  3.1  compares  prediction  accuracy  obtained  by  iterating  equation 
1.1  with  and  without  the  N  term.  Including  N  results  in  prediction 
along  the  Frcnet  circle  while  omitting  its  results  in  constant  accele¬ 
ration  extrapolation.  The  two  predictors  are  compared  for  various  pre¬ 
diction  times  (e.g.,  bullet  times  of  flight)  along  various  coordinated 
turns.  The  table  indicates  that  the  relative  superiority  of  the  Frenet 
circle  predictor  grows  with  target  g's  and  prediction  times.  In  gen¬ 
eral,  it  is  at  least  as  good  as  the  constant  acceleration  extrapolator 
along  noncoordinated  maneuvers. 

Tlicse  ideas  can  be  practically  applied  to  real  systems.  The 
estimator-predictor  algorithms,  equations  1.1  through  1..S,  slightly 
modified  to  include  the  singular  case,  K=0,  have  proven  to  be  complete¬ 
ly  stable  for  circular,  popup  and  dive,  straight  line  and  jinking  tar¬ 
get  ti'aicctories .  The  engineering  expertise  required  to  construct  the 
appropriate  hardware  has  already  been  successfully  demonstrated  in  the 
GLAADS  program.  In  fact,  the  improved  estimator-predictor  could  easily 
be  programmed  into  the  present  GLAADS  computer. 

In  conclusion,  the  plant  noise  mean  and  the  dynamical  plant 
methods  provide  improvements  to  current  state  estimation  technology 
and  offer  an  approach  to  incorj^orating  real  target  constraints  into 
fire  control  logic.  In  terms  of  practical  payoff,  these  methods  will 
enable  the  proposed  Army  Radar  Gun  Air  Defense  System  (ARGAD3)  to  more 
successfully  engage  maneuvering  targets. 


MAGNITUDE  OF  ACCELERATION  ERROR  (M/SEC^ 
SEMI- LOG 


Figure  3.1.  Comparison  of  Estimator  Accuracy. 
Frenet  Circle/Constant  Acceleration 
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